      SUBROUTINE  CDCOEF(C,M,N,IPO,AINT,JWD,D,A,PI,SIGMA,DSUM,FCTR)     00010020
CCC                                                                     00020000
CCC           EIGENVALUES AND EXPANSION COEFFICIENTS OF PROLATE &       00030000
CCC           OBLATE SPHEROIDAL WAVE FUNCTIONS WITH COMPLEX 'C'.        00040000
CCC           IPO = 1  -  PROLATE SYSTEM ;   IPO = -1  -  OBLATE SYSTEM 00050000
CCC           AINT :  INITIAL VALUE OF THE EIGENVALUE.                  00060000
CCC           A  :  FINAL EIGENVALUE CORRECTED BY BOUWKAMP,S METHOD     00070000
CCC           D  :  EXPANSION COEFFUCIENTS - FLAMMER,S NORMALIZATION    00080000
CCC           DSUM  :  SUMMATION OF THE D-COEFFICIENTS                  00090000
CCC           SIGMA  :  NORMALIZATION FACTOR OF THE ANGLE FUNCTIONS     00100000
CCC           MXJ  :    MAXIMUM NUMBER OF THE EXPANSION COEFFICIENTS    00110000
CCC           CALCULATION - COMPLEX DOUBLE PRECISION(COMPLEX*16)        00120000
CCC           C & RESULTS - COMPLEX*16                                  00130000
CCC           EXTERNAL FUNCTION: FF, FACTMM                             00140000
CCC              ORIGINALLY PROGRAMMED BY S.ASANO IN 1972.              00145000   
CCC              MODEFIED IN 1977.                                         
CCC
      PARAMETER   (MXJ=150)                                             00150015
      COMPLEX*16  A,AINT,U1,U2,V1,V2,P1,P2,DELTA,D0,EVEN,ODD,PI,SIGMA,  00160010
     1            DSUM,CONSIG,ALP0,UPLMT,C,C2,A0,DN(MXJ),UN(MXJ),       00170021
     2            CN(MXJ),ALP(MXJ),BET(MXJ),GAM(MXJ),D(MXJ)             00180012
      REAL*8   RM,DM,ENORM,ONORM,R,BETA,FF,FACTMM,FLOAT,Q1,Q2,Q3,Q4     00190021
C                                                                       00200003
      FLOAT(N)= DFLOAT(N)                                               00210003
      BETA(K,I)=FLOAT(I*(I-1))*FLOAT((K+K+I)*(K+K+I-1))/(FLOAT((2*(K+I) 00220000
     +          -3)*(2*(K+I)-1))*FLOAT((2*(K+I)-1)*(2*(K+I)+1)))        00230009
C      
C     CALL ERRSET(207,256,-1,0,0,208)                                   00240000
C                                                                       00250000
          IF(IABS(IPO).EQ.1)  THEN                                      00260002
      C2=C*C*IPO                                                        00270002
      DSUM=(0.D0,0.D0)                                                  00280000
      H=CDABS(C)                                                        00290000
      FCTR= 1.D0                                                        00300019
      IF(H.GE.20.D0.AND.M.GE.17)  FCTR= 1.D+40                          00300122
C                                                                       00300219
      A=AINT                                                            00301019
      L=N-M                                                             00310000
      IFIN0= 40 + 2.7*H + L/5 - M/2                                     00320023
       IF(H.LT.2.0) IFIN0= 20 + 2.0*H + L/5 - M/2                       00320123
       IF(IFIN0.GE.MXJ-1) IFIN0= MXJ-2                                  00330016
      CRIT=1.0E-06                                                      00340000
      IF(H.GE.3.0) CRIT=5.0E-06                                         00360000
      IF(H.GE.15.0) CRIT=1.0E-05                                        00370000
      CONVG= 1.D-10                                                     00380012
      IF(CDABS(A).LE.1.D0) CONVG= 1.D-07                                00390007
CCCC          PARAMETERS - ALPHA,BETA,GAMMA AND RATIO-N                 00400000
      L0=L-2                                                            00410000
        IF(MOD(L,2).EQ.0)  THEN                                         00420003
      INT=2                                                             00430000
      IFIN=2*(IFIN0/2)                                                  00440000
      ALP0=FLOAT(M*(M+1))+0.5D0*C2*(1.D0-FLOAT(4*M*M-1)/FLOAT((2*M-1)*  00450017
     + (2*M+3)))                                                        00460007
        ELSE                                                            00470003
      INT=1                                                             00480003
      IFIN=2*(IFIN0/2)+1                                                00490000
        END IF                                                          00500003
 2002 DO 50 JR=INT,IFIN,2                                               00510000
      RM=DFLOAT(M+JR)                                                   00520000
      DM=DFLOAT(M+M+JR)                                                 00530000
      ALP(JR)=RM*(RM+1.0)+0.5D0*C2*(1.D0- FLOAT(4*M*M-1)/               00540017
     +        ((2.D0*RM-1.D0)*(2.D0*RM+3.D0)))                          00550017
      BET(JR)=C2*C2*BETA(M,JR)                                          00560000
      GAM(JR)=C2*DM*(DM-1.D0)/((2.D0*RM-1.D0)*(2.D0*RM+1.D0))           00570017
       IF(JR.GE.3) THEN                                                 00580003
      UPLMT=BET(JR)/(A-ALP(JR))                                         00590000
      IF(CDABS(UPLMT).LE.CRIT) GO TO 1                                  00600000
       END IF                                                           00610003
   50 CONTINUE                                                          00620003
      JR=IFIN+2                                                         00630000
    1 JW=JR-2                                                           00640000
      CRIT=CRIT/10.                                                     00650000
      IF(JW.LT.10) GO TO 2002                                           00660000
      IF(JW.GE.L) GO TO 300                                             00670000
       GO TO 2002                                                       00680000
CCC           CORRECTION FOR EIGENVALUE  **  BOUWKAMP,S METHOD          00690000
  300 ITERAT=0                                                          00700000
      ISWT=0                                                            00710000
      A0=A                                                              00720000
      LC0=L                                                             00730000
      IS0=2                                                             00740000
      CON=0.8                                                           00750017
      IF(IPO.EQ.-1) CON=1.2                                             00760017
C            ITERATION LOOP                                             00770004
    6 ITERAT=ITERAT+1                                                   00780000
        IF(MOD(L,2).EQ.0)  THEN                                         00790003
      UN(2)=A-ALP0                                                      00800000
      LA=2                                                              00810000
        ELSE                                                            00820003
      UN(3)=A-ALP(1)                                                    00830003
      LA=3                                                              00840000
        END IF                                                          00850003
      JWA=JW-2                                                          00860003
      DN(JW)=BET(JW)/(A-ALP(JW))                                        00870000
      DO 60 IR=LA,JWA,2                                                 00880000
      IRA=JWA+LA-IR                                                     00890000
   60 DN(IRA)=BET(IRA)/(A-ALP(IRA)-DN(IRA+2))                           00900000
      LC=L                                                              00910000
       IF(H.GE.18.0) LC=L+2                                             00920003
       IF(H.GE.38.0) LC=L+4                                             00930003
       IF(CDABS(A0).LT.1.0) LC=L                                        00940003
       IF(L0.LT.1) LC=L+2                                               00950003
       IF(ISWT.EQ.IS0) LC=LC0                                           00960003
      DO 70 KR=LA,LC,2                                                  00970000
   70 UN(KR+2)=A-ALP(KR)-BET(KR)/UN(KR)                                 00980000
      U1=UN(LC)                                                         00990000
      U2=-DN(LC)                                                        01000000
      UPRD=(U1/CDABS(U1))*(U2/CDABS(U2))                                01010000
      IF(UPRD.LE.0.0) GO TO 5050                                        01020000
      ISWT=ISWT+1                                                       01030000
      IF(ISWT.EQ.1.AND.(LC-L).LE.2) LC0=LC+2                            01040000
      IF(ISWT.EQ.1.AND.(LC-L).GE.4) LC0=LC-2                            01050000
      IF(ISWT.EQ.IS0.AND.ITERAT.EQ.IS0) GO TO 5060                      01060000
      GO TO 5050                                                        01070000
 5060 A=A0                                                              01080000
        GO TO 6                                                         01090003
 5050 CONTINUE                                                          01100000
      V2=(0.D0,0.D0)                                                    01110000
      P2=(1.D0,0.D0)                                                    01120000
      DO 90 IQ=LC,JW,2                                                  01130000
      IF(CDABS(P2).LT.1.0D-50) GO TO 90                                 01140000
      P2=P2*DN(IQ)*DN(IQ)/BET(IQ)                                       01150000
      V2=V2+P2                                                          01160000
   90 CONTINUE                                                          01170000
      V1=(1.D0,0.D0)                                                    01180000
      IF(LC.EQ.LA) GO TO 5                                              01190000
      P1=(1.D0,0.D0)                                                    01200000
      LB=LC-2                                                           01210000
      DO 80 KQ=LA,LB,2                                                  01220000
      KK=LA+LB-KQ                                                       01230000
      P1=P1*BET(KK)/UN(KK)/UN(KK)                                       01240000
       IF(CDABS(P1/V1).LT.1.D-50) GO TO 5                               01250006
   80 V1=V1+P1                                                          01260000
    5 DELTA=-(U1+U2)/(V1+V2)                                            01270000
      IF(CDABS(DELTA).LT.1.D-11) GO TO 4                                01290000
      CRITA=CDABS(DELTA/A)                                              01300000
      IF(CRITA.LE.CON) GO TO 502                                        01310000
       GO TO 6                                                          01320004
  502 A=A+DELTA                                                         01330000
        IF(ITERAT.GE.15) GO TO  4                                       01331017
       IF(CRITA.GT.CONVG) GO TO 6                                       01340004
    4 CONTINUE                                                          01350000
C             ITERATION LOOP END                                        01360004
        IF(MOD(L,2).EQ.0)  THEN                                         01370003
      CN(2)=A-ALP0                                                      01380000
      LA=2                                                              01390000
        ELSE                                                            01400003
      CN(3)=A-ALP(1)                                                    01410003
      LA=3                                                              01420000
        END IF                                                          01430003
      JWA=JW-2                                                          01440003
      LB=LC+2                                                           01450000
      CN(JW)=BET(JW)/(A-ALP(JW))                                        01460000
      DO 3020 JA=LB,JWA,2                                               01470000
      JB=JWA+LB-JA                                                      01480000
 3020 CN(JB)=BET(JB)/(A-ALP(JB)-CN(JB+2))                               01490000
       IF(LA.EQ.LC) GO TO 3050                                          01500003
      DO 3030 JC=LA,LC,2                                                01510000
 3030 CN(JC+2)=A-ALP(JC)-BET(JC)/CN(JC)                                 01520000
 3050 CONTINUE                                                          01530000
CCC           EXPANSION COEFFICIENTS - D                                01540000
      M2=M+M                                                            01550000
        IF(MOD(L,2).EQ.0) THEN                                          01560003
      LD=4                                                              01570000
      D(2)=CN(2)/GAM(2)                                                 01580000
      R=-1.0D0                                                          01590000
       IF(MOD(L/2,2).EQ.0) R=1.0D0                                      01600003
      ENORM=R/FF(M,N)                                                   01610000
      EVEN= FACTMM(M)                                                   01620018
        ELSE                                                            01630003
      LD=5                                                              01640003
      D(3)=CN(3)/GAM(3)                                                 01650000
      LR=(L-1)/2                                                        01660000
      R=-1.D0                                                           01670000
       IF(MOD(LR,2).EQ.0) R=1.D0                                        01680003
      ONORM=R/FF(M,N)                                                   01690000
      ODD= 0.5D0*FACTMM(M+1)                                            01700018
        END IF                                                          01710003
      LE=LD-2                                                           01720003
      DO 15 I=LD,JW,2                                                   01730000
      IN=I                                                              01740000
       IF(CDABS(D(I-2)).LT.1.D-74) GO TO 13                             01750022
   15 D(I)=CN(I)/GAM(I)*D(I-2)                                          01760000
      JWD=JW                                                            01770000
      GO TO 14                                                          01780000
   13 JWD=IN                                                            01790000
   14 DO 17 I=LE,JWD,2                                                  01800000
       IRM=I+M                                                          01801022
       CALL  SUBFF(M,IRM,FFD,C0FF)                                      01802022
        IF(MOD(L,2).EQ.0)  THEN                                         01810003
      IA=I/2                                                            01820000
      R=1.D0                                                            01833022
       IF(MOD(IA,2).EQ.1) R=-1.D0                                       01840003
      EVEN=EVEN + C0FF*R*D(I)/FFD                                       01860022
        ELSE                                                            01870003
      IB=(I-1)/2                                                        01880003
      R=1.D0                                                            01890000
       IF(MOD(IB,2).EQ.1) R=-1.D0                                       01900003
      ODD=ODD + C0FF*R*D(I)/FFD                                         01920022
        END IF                                                          01930003
   17 CONTINUE                                                          01940003
         IF(MOD(L,2).EQ.0)  THEN                                        01950003
      D0=ENORM/EVEN                                                     01960000
      D(1)=D0                                                           01970000
      DO 25 J=2,JWD,2                                                   01980000
   25 D(J+1)=D0*D(J)                                                    01990000
         ELSE                                                           02000003
      D(1)=ONORM/ODD                                                    02010003
      DO 35 J=3,JWD,2                                                   02020000
   35 D(J)=D(1)*D(J)                                                    02030000
         END IF                                                         02040003
CC         SUM OF D-COEF. , NORMALIZATION CONST. AND FACTOR PI          02050007
      DSUM=(0.D0,0.D0)                                                  02060003
      PI=(0.D0,0.D0)                                                    02070000
       JA=JWD                                                           02080007
       IF(MOD(L,2).EQ.0) JA=JWD+1                                       02090003
      DO 21 JB=1,JA,2                                                   02100000
      JJ=JA+1-JB                                                        02110000
       J=JJ-1                                                           02120007
       IF(MOD(L,2).EQ.1) J=JJ                                           02130003
      DSUM=DSUM+D(JJ)                                                   02140000
      Q1= 1.D0                                                          02150021
      Q2= 1.D0                                                          02151021
        IF(M.GE.1)  THEN                                                02160017
      DO 19 LL=1, M                                                     02170017
      Q1= Q1*DFLOAT(J+LL)                                               02180017
   19 Q2= Q2*DFLOAT(J+LL+M)                                             02181017
        END IF                                                          02190011
      PI=PI+ D(JJ)*Q1/FCTR*Q2                                           02200021
   21 CONTINUE                                                          02210000
      SIGMA=(0.D0,0.D0)                                                 02220000
      DO 31 JJ=1,JWD,2                                                  02230000
       J=JJ-1                                                           02240007
       IF(MOD(L,2).EQ.1) J=JJ                                           02250003
      Q3= 1.D0                                                          02260021
      Q4= 1.D0                                                          02261021
        IF(M.GE.1)  THEN                                                02270017
      DO 38 II=1, M                                                     02280017
      Q3= Q3*DFLOAT(J+II)                                               02290017
   38 Q4= Q4*DFLOAT(J+II+M)                                             02291017
        END IF                                                          02300007
      CONSIG= 2.D0/DFLOAT(M2+J+J+1)*D(JJ)*Q3/FCTR*D(JJ)*Q4              02320021
      SIGMA=SIGMA+CONSIG                                                02330000
       IF(J.GT.JWD/2.AND.CDABS(CONSIG/SIGMA).LT.1.D-25) GO TO 33        02340019
   31 CONTINUE                                                          02350000
   33 CONTINUE                                                          02360000
        RETURN                                                          02370002
C* 30 WRITE(6,1080)                                                     02380017
C       RETURN                                                          02390017
C1080 FORMAT(1H0,20X,'*  CDCOEF  -  RETURN WITHOUT CALCULATION  *  EIGEN02400017
C*   2VALUE IS NOT DETERMINED ' )                                       02410017
          ELSE                                                          02420002
      WRITE(6,1050)                                                     02430002
        RETURN                                                          02440002
 1050 FORMAT(1H0,20X,'*  CDCOEF  -  RETURN WITHOUT CALCULATION  *  IPO  02450002
     1SHOULD BE 1 OR -1' )                                              02460002
          END IF                                                        02470002
      END                                                               02480000
